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SUMMARY 


A method is presented for using the results of optimal control theory to 
synthesize a feedback filter. The feedback filter is used to force the output 
of the filtered frequency response to match that of a desired optimal frequency 
response over a finite frequency range. This matching is accomplished by 
employing a nonlinear programing algorithm to search for the coefficients of 
the feedback filter that minimize the error between the optimal frequency 
response and the filtered frequency response. The method is applied to the 
synthesis of an active flutter-suppression control law for an aeroelastic wind- 
tunnel model. It is shown that the resulting control law suppresses flutter 
over a wide range of subsonic Mach numbers. The study indicates that this is 
a promising method for synthesizing practical control laws using the results of 
optimal control theory. 


INTRODUCTION 

In modern control theory, a number of methods have been developed to 
design control systems which require multiple variables (so-called state 
variables) to describe the state or condition of the system. Optimal regu- 
lator theory (ref. 1) is the most widely used method for determining optimal 
control laws that feed back all of the state variables (full-state feedback) . 
Most applications of optimal regulator theory have been limited to problems 
with only a few state variables. Recently, however, this theory has been 
applied to aeroelastic problems that are generally characterized by a large 
number of state variables (ref. 2) . 

Several researchers have applied optimal regulator theory to active flut- 
ter suppression. The two-dimensional flutter problem was investigated by 
Lyons et al. (ref. 3) and Edwards (ref. 4) . Their work was extended in refer- 
ence 5 by considering a more complete mathematical model incorporating three- 
dimensional unsteady aerodynamics. It is shown in reference 5 that optimal 
regulator theory provides full-state feedback control laws that are very 
attractive for flutter suppression. 

A major problem in the application of optimal regulator theory is the 
conversion of the theoretical (full-state feedback) control law into a prac- 
tical control law (that is, one that can be readily implemented). This problem 
arises because optimal regulator theory normally requires that all state vari- 
ables be available for feedback. Direct measurement of all feedback states for 
an aircraft wing, for example, is not possible, and a method is required to 
obtain a practical control law from the limited sensor measurements available. 
Konar et al. (ref. 6) developed such a method by adjusting the optimal full- 
state feedback gains to be compatible with sensor measurements. This method 
uses a numerical search algorithm in the time domain that adjusts feedback 
gains while maintaining the least amount of increase in the quadratic objective 


function. This method has been applied to studies involving load alleviation 
for the C-5A (ref. 2) . 

The purpose of this paper is to describe a different approach that uses 
frequency domain techniques to obtain a practical control law. This approach 
employs a transfer- function matching technique developed by Coffey (ref. 7) . 
This method uses a gradient optimization algorithm to design a feedback filter 
that forces the open-loop frequency response of the system to match a desired 
open-loop frequency response. In reference 7, it is assumed that the desired 
open-loop frequency response is known. The present approach defines the 
desired open-loop frequency response as that of the optimal full-state feedback 
system. Also given is a brief description of the method employed to develop 
the equations in terms of state variables. The technique is applied to the 
synthesis of a practical control law for active flutter suppression of an aero- 
elastic wind-tunnel model. 
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SYMBOLS 

ith denominator coefficient of actuator transfer function 

numerator coefficient of actuator transfer function 

reference length, m 

steamwise local chord, m 

error function in frequency domain 

frequency, Hz 

factor in feedback filter 

gravitational constant, 9.80 m/sec^ 

feedback-filter transfer function 

quadratic optimization cost function 
scalar gain for feedback filter 

reduced frequency, — 

2U 

number of frequencies 

number of aerodynamic lag terms 

flutter dynamic pressure, kPa 
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<*o 


free- stream dynamic pressure. 



kPa 


s Laplace variable 

t time, sec 

free- stream velocity, m/sec 
u control input 

u practical control input 

u optimal control input 

Wg vertical gust velocity, m/sec 

y output variable 

z vertical displacement, m 

3 m aerodynamic lag terms 

5 control surface position, positive down, deg 

5 C actuator command, deg 

£ damping ratio 

w circular frequency, rad/sec 

to n natural frequency, rad/sec 

Matrices: 


[a] 

[A a ] 

[A^ 

[A v ] 

[ A(ik) ] 

[b] 

[B a ] 

[b v ] 


total-system dynamics matrix 
actuator dynamics matrix 
real aerodynamic matrix coefficients 
vehicle dynamics matrix 

generalized aerodynamic-force matrix due to wing motion 
total-system control distribution matrix 
actuator control distribution matrix 


vehicle control distribution matrix 



[C] 

tc a ] 

[c s ] 

[c v ] 

Ic^j 

[D v ] 

[D (ik) ] 
{G (ik) } 
[I] 

[K] 

[K s ] 

[Mg] 

[P] 

[Q] 

{q> 

[R] 

{u} 

{u a } 

{u v } 

{u} 

{X} 

{x a } 

{X v } 

(Y> 

{Y a } 

{Y v } 


total-system state-coefficient output matrix 
actuator state-coefficient output matrix 
generalized damping matrix 
vehicle state-coefficient output matrix 
row matrix of mode-shape amplitudes 
vehicle input-coefficient output matrix 

generalized aerodynamic-force matrix due to control surface rotation 

generalized aerodynamic-force vector due to vertical gust velocity 

identity matrix 

optimal gain matrix 

generalized stiffness matrix 

generalized mass matrix 

Riccati matrix 

output weighting matrix 

generalized coordinate vector 

control weighting matrix 

total- system input vector 

actuator input vector 

vehicle input vector 

optimal control input vector 

total-system state vector 

actuator state vector 

vehicle state vector 

total- system output vector 

actuator output vector 

vehicle output vector 
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{ 6 } 


control surface position vector 


Subscr ipts: 


f 

flutter 


I 

imaginary 

part of complex value 

R 

real part 

of complex value 

rms 

root-mean- 

-square value 


Dots over symbols denote derivatives with respect to time. 

EQUATIONS OF MOTION 

The equations of motion for a flexible vehicle may be expressed in matrix 
form as 

w g 

[ M S I {q} + [c s ]{q} + [k s 1 { q} + q^t A( ik) ] { q) + q^f D(ik)] {6> = q^Gfik)} — (1) 

where 

[m s ] generalized mass matrix 

[Cg] generalized damping matrix 

[k s ] generalized stiffness matrix 

[ A( ik) ] generalized aerodynamic-force matrix due to wing motion 

[D(ik)] generalized aerodynamic-force matrix due to control surface rotation 

{G(ik)} generalized aerodynamic- force vector due to vertical gust velocity 

In flutter analyses where unsteady aerodynamics are of major importance, 
the aerodynamic-force matrices are normally represented as tabular functions 
of reduced frequency k = cay^U. with the unsteady aerodynamic forces in this 
form, the equations of motion cannot be written as a set of first-order dif- 
ferential (state-space) equations. The equations can be cast into state-space 
form, however, by use of aerodynamic approximation functions (refs. 8 and 9). 

Aerodynamic Approximation Functions and State-Space Equations 

The techniques for developing the aerodynamic approximation equations and 
the state-space equations are similar to those described in references 5 and 10. 
The aerodynamic approximation technique involves fitting the curve of a mathe- 
matical function of ik to the aerodynamic forces. Each element of the 
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aerodynamic-force matrix is fit with an approximation function of the form 


„ M [A m+2 ] (ik) 

[ A(ik)] = [Aq] + [A^ (ik) + [A 2 ] (ik ) 2 + 2lI ~ 7 (2) 

m=l &m + ( lk > 

The aerodynamic lag terms $ m are somewhat arbitrary but are usually chosen 
to lie somewhere within the range of interest of the reduced frequencies. The 
numerator matrices are computed to give a least-square error for the values of 
k at which the aerodynamic forces are known. 

As described in reference 5, substituting equation (2) into equation (1) 
and equating derivatives to the powers of the Laplace operator 
(s = io) = (2U/c)ik) yields the equations of motion for the basic vehicle. 
Written in standard state-space form, 

{x v } = [a v Hx v } + [b v Hu v m 

> (3) 

{y v } = [c v Hx v }+ [d v Hu v }J 

The elements of the vehicle state vector {x v } are the generalized coordinates 
and their derivatives. The vehicle input vector {u v } consists of the control 
surface displacement, the gust disturbance, and their derivatives. The order 
of the derivatives depends on the number of lag terms used in the aerodynamic 
approximation. 


Actuator Models 

To include the effect of actuator dynamics during design of the optimal 
control law, the actuator model is described in state-space form and then 
interconnected to the basic vehicle equations. Actuator models are generally 
represented by transfer functions. Consider an actuator transfer function of 
the following form: 

<S(s) b o 

7 = (4) 

M s > s" + a n _-] s n_1 + . . . a 0 

By cross multiplying and equating derivatives to the powers of s, an nth order 
differential equation of the actuator can be written 


d5 

+ a n-l 

dt n 


d<5 

+ 

dt n_1 


a 0 6 = Me 


(5) 


By making the substitutions that X^ = d* -1 6/dt*" 1 , a set of n first-order 
differential equations can be written. Thus, 
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The output equations of the actuator are written in the form 

(Y a > = [C a Hx a } (7) 


where 


[C a l 


1 0 0 
0 1 0 


The interconnection of the actuator to the basic vehicle involves equating 
the output of the actuator to the input of the basic vehicle, as illustrated in 
the following sketch: 


Total System 


{u a } 



{Y a } = {u v } 



(Y v > 
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By using the relationship {u v } = {Y a }, the total-system equations can be 
written as 



( 8 ) 


(9) 


SYNTHESIS OF AN ACTIVE FLUTTER-SUPPRESSION CONTROL LAW 


Regulator Theory 


Optimal regulator theory provides for the minimization of a quadratic cost 
function of the output and control vectors (ref. 1). To find the optimal full- 
state feedback control law, the quadratic cost function 


J 


^ ^{y} T [qHy>+ {u } T [r] (u ^ dt 


( 10 ) 


is minimized. This leads to the optimal control law 


(u> 

Ck 1 


= {u} = [kHx} 

= - [r] - 1 [ B ] T [p] 


( 11 ) 


For the time invariant (constant-coefficient differential equations) condition, 
[p] is the steady state solution of the matrix Riccati equation. 

[P] [A] + [A] T [P] - [P] [B] [R] _1 [ B ] T [P] + [Q] * 0 (12) 

The application of regulator theory is an iterative process of selecting the 
appropriate cost function through changes in the weighting matrices [q] and 
[r] . The procedure can be summarized as follows: 
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Step 1 : Define the output ({y}) and input ({u}) vectors that relate 

to the performance goals. (For example, minimum control 
surface displacement . ) 

Step 2: Select initial weighting matrices [q] and Dr ] . (For 

example, see refs. 1 to 5.) 

Step 3: Solve equations (11) and (12) for the optimal gains (fel), 

thus minimizing the quadratic cost function J. 

Step 4: Evaluate the design, and adjust weighting matrices until 

performance goals are met. 

The optimal control law requires the capability to feed back all of the 
state variables. Since the description of the system involves modal or gener- 
alized coordinates in addition to physical coordinates, direct measurement of 
all state variables is not feasible. Only a linear combination of the state 
variables can be measured. Therefore, a method is employed which uses the 
available measurements. 


Process for Design of Practical Control Law 

The process described herein for the design of a practical control law is 
performed in the frequency domain and attempts to match a desired open-loop 
frequency response. The desired open-loop frequency response is that of the 
optimal full-state feedback system. The design process involves finding the 
coefficients of a feedback filter H(s) that minimize the deviation of the 
open-loop frequency response (u/u) (ioo) from the optimal open-loop frequency 
response (u/u)(ioo). (See fig. 1.) Figure 2 is a Nyquist diagram illustrat- 
ing the results of the design process for a practical control law. The objec- 
tive is to make the deviation from the optimal system small. If this objective 
is met, the performance of the practical system will be similar to that of the 
optimal system. 

Error function .- The error function is defined as the difference of 

(u/u) (ico) and (u/u) (iw) over a set of frequency points 

(i = 1, 2, . . L) , for which a close fit is desired. The closeness of 

(u/u) (iw) to (u/u) (ico) can be described mathematically by 
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-1 
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u 

-(iw)* - 

-(ia))* 


-(ia» 

- -(ico) 

u 

u 


u 

u 


E = 


(13) 


where the asterisk (*) denotes the complex conjugate, with 



( 14 ) 


Feedback-filter design variables .- The form of the feedback filter 
to be used is 


m 

n (s 2 + 2CiOJ ni s + co 2 i) 
i=l 

H(s) = K f f(s) (15) 

n 

n (s 2 + 2CjW n jS + 0)2j) 
j=l 

The design variables are the gain K.f, the damping ratios C, and the frequen- 
cies a) n (of each second-order factor in eq. (15)). The function f(s) is 
included to help achieve any desired characteristics of the feedback filter, 
such as high-frequency roll-off. The procedure can be described as follows: 

Step 1 : Compute the optimal frequency response u/u. 

Step 2: Compute the frequency response between the output y 

and the control u. 

Step 3: Choose the initial number of numerator factors m and 

denominator factors n of the filter. 

Step 4: Choose a f (s) to incorporate any desired characteristics 

of the filter, such as high-frequency roll-off. 

Step 5: Minimize the error function E by using an optimization 

algorithm such as that of Davidon (ref. 11) and Fletcher 
and Powell (ref. 1 2) . 

Step 6: Examine the practical open-loop frequency response (Nyquist 

diagram) to establish any possible changes to f (s) . 

Step 7: If any changes to f(s) are established, repeat step 5. 

Step 8: Repeat steps 3 through 7 for a family of m and n. 
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Step 9: Select the values of m, n, and f(s) that provide the 

smallest value of E. 


APPLICATION OF TECHNIQUE 

The methodology described previously is applied to the synthesis of an 
active flutter-suppression control law for an aeroelastic wind-tunnel model. 
The model geometry is shown in figure 3. The model consists of a cantilever 
wing with a 20-percent-chord, trailing-edge control surface located between 
the 76-percent and 89-percent semispan stations. The flutter-suppression 
sensor (accelerometer) is located at the 60-percent-chord and 92-percent- 
semispan station. The first eight elastic modes 1 are used as generalized 
coordinates covering a frequency range from 5.23 Hz to 118.15 Hz. The calcu- 
lations necessary to determine the coefficients of the equations of motion 
(eq. (1)) are described in detail in reference 10. 


Basic Wing Characteristics 

Each of the aerodynamic terms is approximated in the s-plane through the 
use of equation (2) with M = 2. The 8 m terms are varied until an acceptable 
curve fit is found. This resulted in the B m terms being selected as 0.225 
and 0.500 for all aerodynamic terms. In figure 4, one of the calculated 
oscillatory aerodynamic terms is compared with the approximation function at 
Mach =0.9. In general, all of the aerodynamic terms have a good curve fit. 

To validate the mathematical model further, the flutter boundary of the 
model without the flutter-suppression system (FSS off) is calculated. For a 
specific Mach number, the characteristic roots of equation (3) are found for a 
series of dynamic pressures. The dynamic pressure at which the real part of 
one of the roots becomes zero is the flutter dynamic pressure. Shown in fig- 
ure 5 are the dynamic-pressure root loci at Mach = 0.9. Calculations of the 
dynamic-pressure root loci were also performed at Mach = 0.6, 0.7, and 0.8 
to establish the FSS-off flutter boundary shown in figure 6. A comparison of 
the experimental results reported in reference 10 with these analyses indicates 
good agreement. 

In addition to verifying the aerodynamic approximations with respect to 
the basic flutter characteristics, the transfer function between acceleration 
and control surface deflection (FSS off) is compared to that using the origi- 
nal, oscillatory aerodynamic forces. In this manner, the approximations for 
the control surface aerodynamic terms are verified. The gain and phase curves 
for this transfer function using the original oscillatory aerodynamics (desig- 
nated as k-plane) and using the aerodynamic approximation function (designated 
as s-plane) are shown in figure 7 for comparison. Good agreement is indicated 
in both gain and phase. 


^The two inplane modes (3 and 8) of reference 10 are omitted in this 
study. 
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Optimal and Practical Control Laws 


The design condition for the FSS is chosen to be a 44-percent increase in 
flutter dynamic pressure at Mach = 0.9. The FSS is required to increase the 
flutter dynamic pressure by at least 44 percent at each of four Mach numbers 
(0.6, 0.7, 0.8, and 0.9) and, at the 44-percent margins, to exhibit ±6 dB gain 
margins and ±30° phase margins. In addition, the control surface activity can 
not exceed 6° (rms) and 600°/sec (rms). The optimal control laws are deter- 
mined by solving equations (1 1 ) and (1 2) and thus satisfying equation (1 0) . 
Equation (12) is solved using the computational algorithm developed by Vaughan 
(ref. 13) as coded in reference 14. For the quadratic optimization, zero-state 
weighting ([q] = 0) is selected since this yields a set of gains that are 
"cheapest" (ref. 15) in terms of control input amplitude. This optimal full- 
state feedback control law leaves all stable eigenvalues unchanged and relo- 
cates the unstable eigenvalues to their mirror image in the left half plane. 
Once the full-state feedback gain matrix is determined, the optimal Nyquist 
diagram for the single input system is constructed by solving 


LkJ L itoi - a] Ib> (16) 

u(iw) 

Equation (16) is solved for a series of frequencies from 1 rad/sec to 
301 rad/sec at increments of 3 rad/sec. The resulting Nyquist diagram 
(fig. 8) is a counterclockwise circle of radius unity centered on the (-1 , 0) 
point. The full-state feedback control law provides gain and phase margins of 
-6 dB and ±60°, respectively. Note that the Nyquist diagram crosses the real 
axis only once (-6 dB) and, therefore, has infinite positive gain margin. 

To synthesize a practical control law from the full-state feedback control 
law, the frequency response between sensor output and control surface input is 
established first. Using the accelerometer location shown in figure 3, the 
output frequency response is calculated by 


Z(1W) , , r , -l r ! 

= lC<hJ tiWI - Al {B> (17) 

u(iw) 

where lc<f,J is a row matrix of mode-shape amplitudes at the sensor location. 
The problem is to find the coefficients of the feedback filter H(iiO) that 
satisfy 


z(iw) u ( i w) 

x H(iW) = (18) 

u(iw) u(iw) 


over a finite set of frequency points (i = 1, . . ., L). In the frequency 

plane, equation (15) has the form 
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H(i(d) = K f 


n (ico) 2 + 2 CiW ni (iu) + w&jl 

i=l u 

n (i 0 )) 2 + 2 Cj(*) n j (iw) + Wnj 

j=l L 


(19) 


Employing a trial and error approach, various values of m and n are tried, 
and the combination that provides the smallest value of E (m = 2, n = 2) is 
determined. The low-frequency (3 rad/sec to 60 rad/sec) portion of the open- 
loop frequency response (not shown) indicated the need for an integrator. 
Therefore the f(io)) is selected to be 1 / ( ia) ) . The error function is again 
minimized with the previously determined values of m and n and the l/(iw) 
in the filter. Figure 9 is the Nyquist diagram resulting from the minimization 
process. Since the optimization algorithm did not result in a -6 dB gain mar- 
gin, the gain Kf was increased until a gain margin of -6 dB was achieved 
which resulted in the following control law: 


6 2214 s 2 + 2(0.127) (121 .21 )s + (121. 21 ) 2 

z s S 2 + 2(0.962) (297. 62)s + (297. 62) 2 

s 2 + 2(0.088) (269. 14)s + (269 . 1 4 ) 2 /deg\ 
x (20) 

s 2 + 2(0.964) (294. 91)s + (294.91) 2 \9 / 

Figure 10 is a Nyquist diagram obtained by using the control law defined by 
equation (20) . 


Results 

The performance of the control law in terms of increased flutter dynamic 
pressure, root-mean- square (rms) values of the control activity in turbulence, 
and gain/ phase margins are examined at Mach = 0.6, 0.7, 0.8, and 0.9. All cal- 
culations are performed using a modified version of the computer program 
described in reference 16. A summary of these results is presented in table 1. 

Flutter characteristics .- Dynamic-pressure root locus calculations are 
performed to establish the FSS-on flutter boundary. Shown in figure 11 are 
the FSS-on root loci at Mach = 0.9. The control law increases the damping 
of the flutter mode while having very little effect on the other modes. The 
increased damping delays flutter onset until q^ = 9.863 kPa, which is a 
96-percent increase in flutter dynamic pressure. Large increases in flutter 
dynamic pressure are predicted at the other Mach numbers as illustrated by the 
FSS-on flutter boundary presented in figure 12. 

Control surface activity .- Control surface activity in turbulence is 
determined using power-spectral-density (PSD) analyses similar to that 
described in reference 10. A Von Karman gust spectrum with a characteristic 
length of 30.48 m is used to simulate turbulence distribution and intensity in 
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the wind tunnel. The variations of <S rms and 6 rms with dynamic pressure at 
Mach = 0.6, 0.7, 0.8, and 0.9 are shown in figure 13. The largest values of 

<5 rms (5.5°) and 6 rms (270°/sec) occur at the maximum dynamic pressure inves- 

tigated (q^ = 10.773 kPa at Mach = 0.6). As indicated, 6 rms and S rms are 
only slightly dependent on Mach number but highly dependent on dynamic pressure. 

Gain and phase margins .- Nyquist diagrams are constructed at all four Mach 
numbers to establish gain and phase margins. The gain margins at the Mach 0.9 
design condition are -6.27 and +13.60 dB with phase margins of -58.8° and 41.0°. 
(See fig. 10.) The Nyquist diagrams at the other Mach numbers (not shown) are 
similar in character to that at Mach 0.9. The gain and phase margins at all 
four Mach numbers are presented in table 1 for a dynamic pressure 44 percent 
above the FSS-off flutter boundary. Gain margins of approximately ±6 dB are 
exhibited at all four Mach numbers; however, the positive phase margins at the 
lower Mach numbers are less than the +30° requirement. 


CONCLUDING REMARKS 

A method is presented for synthesizing feedback control laws using optimal 
regulator theory. The method is applied to the synthesis of an active flutter- 
suppression control law for an aeroelastic wind-tunnel model. Important results 
of the study are: 

1. It is shown that a practical flutter-suppression control law can be 
synthesized by use of a gradient optimization algorithm to design a feedback 
filter which minimizes the difference between the filtered frequency response 
and the optimal frequency response. 

2. Application of the method to a wind-tunnel model provides a control law 
which is shown by analysis to be capable of increasing flutter dynamic pressure 
by at least 44 percent over a range of Mach numbers from 0.6 to 0.9. 


Langley Research Center 

National Aeronautics and Space Administration 
Hampton, VA 23665 
June 4, 1979 
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TABLE 1 


Mach 

number 

<3f 

(FSS off) , 
kPa 

<If 

(FSS on) , 
kPa 

(percent increase] 

0.60 

7.305 

1 2.832 
(75) 

0.70 

6.655 

i 

11.730 

(78) 

0.80 

5.906 

1 0.294 



(76) 

0.90 

5.027 

9.863 



(96) 


*Design point. 


LAW PERFORMANCE 


1 .44q , 

. 00 

kPa 

^rms' 

deg 

• 

6rms' 

deg/sec 

Gain 

margin, 

dB 

Phase 

margin 

deg 

10.773 

5.5 

270 

-5.89 

+9.90 

-33.0 

+16.8 

9.815 

5.1 

258 

-6.24 
+1 0.75 

-32.4 

+20.2 

8.618 

4.5 

238 

-6.49 
+1 2.00 

-41 .0 
+26.6 

*7.661 

4.0 

209 

-6.27 
+1 3.60 

-58.8 
+41 .2 


(a) Optimal. (b) Practical 

Figure 1 Block diagrams of optimal and practical control laws. 
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Figure 3.- Model geometry. 
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Figure 5.- Dynamic-pressure root loci; FSS off, Mach 
(Arrows indicate increasing dynamic pressure.) 
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Figure 7.- Bode diagram of transfer function between acceleration and control surface input 
(FSS off); q_ = 7.661 kPa, Mach = 0.9. (Vertical axis is in decibels.) 
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Figure 11.- Dynamic-pressure root loci; FSS on, Mach = O.S 
(Arrows indicate increasing dynamic pressure.) 
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